library("rstudioapi")     
setwd(dirname(getActiveDocumentContext()$path))
#setwd("../")
getwd()
library("readxl") 
require(miceadds)
library(data.table)
library(tidyverse)
library(DescTools)
library(htmlTable)
library(stargazer)
library(estimatr)
library(tidyverse)
library(geosphere)
library(DescTools)
require(dplyr)

load("FinalFinal20062023Full1-6.RData")
load("FinalFinalSP20062023Full1-6.RData")

# read data:
#----- Calculations 
combo = final_finalV2

#Total number of subjects
subj_total <- nrow(combo)

#Total number of subjects per treatment
subj_treatment <- combo %>%
  group_by(individual_treatment) %>%
  summarise(n = n())

subj_t1 <- subj_treatment$n[subj_treatment$individual_treatment == "Placebo"][1]
subj_t2 <- subj_treatment$n[subj_treatment$individual_treatment == "CDC Health"][1]
subj_t3 <- subj_treatment$n[subj_treatment$individual_treatment == "Low Cash"][1]
subj_t4 <- subj_treatment$n[subj_treatment$individual_treatment == "High Cash"][1]

# show numbers: 
subj_total; subj_t1; subj_t2; subj_t3; subj_t4
row1 = (c(subj_t1, subj_t2, subj_t3, subj_t4,subj_total))

# R=1 
combo$R = ifelse(is.na(combo$p_ii) & is.na(combo$p_iii), 0, 
                 ifelse(combo$p_ii == 1 |  combo$p_iii == 1, 1,0))

comboR = combo[which(combo$R == 1),]
comboR0 = combo[which(combo$R == 0),]


#Total number of subjects
subj_totalR <- nrow(comboR)

#Total number of subjects per treatment
subj_treatmentR <- comboR %>%
  group_by(individual_treatment) %>%
  summarise(n = n())

subj_t1R <- subj_treatmentR$n[subj_treatmentR$individual_treatment == "Placebo"][1]
subj_t2R <- subj_treatmentR$n[subj_treatmentR$individual_treatment == "CDC Health"][1]
subj_t3R <- subj_treatmentR$n[subj_treatmentR$individual_treatment == "Low Cash"][1]
subj_t4R <- subj_treatmentR$n[subj_treatmentR$individual_treatment == "High Cash"][1]

subj_totalR; subj_t1R; subj_t2R; subj_t3R; subj_t4R
row2 = (c(as.numeric(subj_t1R), as.numeric(subj_t2R), as.numeric(subj_t3R), as.numeric(subj_t4R), as.numeric(subj_totalR)))

# Summary Table 1 generated:
sum_tab1 = as.data.frame(cbind((row1), (row2)))
colnames(sum_tab1) = c("Phase 1 - Baseline", "Phase 2 & 3 - Follow-up")
sum_tab1$Change = as.numeric(sum_tab1[,2]- sum_tab1[,1])
sum_tab1$"Percent Change" = as.numeric(round(sum_tab1$Change / sum_tab1[,1] * 100, 2))
sum_tab1 # quick view of the table
cols1 = colnames(sum_tab1)
sum_tab1 <- data.frame(lapply(sum_tab1, as.character), stringsAsFactors=FALSE)
colnames(sum_tab1) = cols1
rownames(sum_tab1) = c("Placebo", "Health", "Low Cash", "High Cash", "Total Sample")

require(stargazer)
stargazer(sum_tab1,type = "latex", summary = FALSE, out="ExtData_Tables3_PanelA.tex")


# Summary Table generated: 
# with vaccine intention: 

# vaccine intention phase 1 combo 
#Vaccine initial intentions total
v_init_int_total <- combo %>%
  group_by(vaccine_intention) %>%
  summarise(n = n())

vii_total_1 <- v_init_int_total$n[v_init_int_total$vaccine_intention == 1][1]/sum(v_init_int_total[which(v_init_int_total$vaccine_intention %in% c(0,1)),]$n)*100

#Vaccine initial intentions per treatment
v_init_int_treatment <- combo %>%
  group_by(individual_treatment, vaccine_intention) %>%
  summarise(n = n())

vii_pla_1 <- v_init_int_treatment$n[v_init_int_treatment$vaccine_intention == 1 & v_init_int_treatment$individual_treatment == "Placebo"][1]/sum(v_init_int_treatment[which(v_init_int_treatment$individual_treatment=="Placebo"&v_init_int_treatment$vaccine_intention %in% c(0,1)),]$n)*100
vii_cdc_1 <- v_init_int_treatment$n[v_init_int_treatment$vaccine_intention == 1 & v_init_int_treatment$individual_treatment == "CDC Health"][1]/sum(v_init_int_treatment[which(v_init_int_treatment$individual_treatment=="CDC Health"&v_init_int_treatment$vaccine_intention %in% c(0,1)),]$n)*100
vii_lc_1 <- v_init_int_treatment$n[v_init_int_treatment$vaccine_intention == 1 & v_init_int_treatment$individual_treatment == "Low Cash"][1]/sum(v_init_int_treatment[which(v_init_int_treatment$individual_treatment=="Low Cash"&v_init_int_treatment$vaccine_intention %in% c(0,1)),]$n)*100
vii_hc_1 <- v_init_int_treatment$n[v_init_int_treatment$vaccine_intention == 1 & v_init_int_treatment$individual_treatment == "High Cash"][1]/sum(v_init_int_treatment[which(v_init_int_treatment$individual_treatment=="High Cash"&v_init_int_treatment$vaccine_intention %in% c(0,1)),]$n)*100

row3 = cbind(vii_pla_1, vii_cdc_1, vii_lc_1, vii_hc_1, vii_total_1)


# vaccine intention R1
v_init_int_totalR <- comboR %>%
  group_by(vaccine_intention) %>%
  summarise(n = n())

vii_total_1R <- v_init_int_totalR$n[v_init_int_totalR$vaccine_intention == 1][1]/sum(v_init_int_totalR[which(v_init_int_totalR$vaccine_intention %in% c(0,1)),]$n)*100

#Vaccine initial intentions per treatment
v_init_int_treatmentR <- comboR %>%
  group_by(individual_treatment, vaccine_intention) %>%
  summarise(n = n())

vii_pla_1R <- v_init_int_treatmentR$n[v_init_int_treatmentR$vaccine_intention == 1 & v_init_int_treatmentR$individual_treatment == "Placebo"][1]/sum(v_init_int_treatmentR[which(v_init_int_treatmentR$individual_treatment=="Placebo"&v_init_int_treatmentR$vaccine_intention %in% c(0,1)),]$n)*100
vii_cdc_1R <- v_init_int_treatmentR$n[v_init_int_treatmentR$vaccine_intention == 1 & v_init_int_treatmentR$individual_treatment == "CDC Health"][1]/sum(v_init_int_treatmentR[which(v_init_int_treatmentR$individual_treatment=="CDC Health"&v_init_int_treatmentR$vaccine_intention %in% c(0,1)),]$n)*100
vii_lc_1R <- v_init_int_treatmentR$n[v_init_int_treatmentR$vaccine_intention == 1 & v_init_int_treatmentR$individual_treatment == "Low Cash"][1]/sum(v_init_int_treatmentR[which(v_init_int_treatmentR$individual_treatment=="Low Cash"&v_init_int_treatmentR$vaccine_intention %in% c(0,1)),]$n)*100
vii_hc_1R <- v_init_int_treatmentR$n[v_init_int_treatmentR$vaccine_intention == 1 & v_init_int_treatmentR$individual_treatment == "High Cash"][1]/sum(v_init_int_treatmentR[which(v_init_int_treatmentR$individual_treatment=="High Cash"&v_init_int_treatmentR$vaccine_intention %in% c(0,1)),]$n)*100

row4 = cbind(vii_pla_1R, vii_cdc_1R, vii_lc_1R, vii_hc_1R, vii_total_1R)
row4

# vaccine intention R0
v_init_int_totalR0 <- comboR0 %>%
  group_by(vaccine_intention) %>%
  summarise(n = n())

vii_total_1R0 <- v_init_int_totalR0$n[v_init_int_totalR0$vaccine_intention == 1][1]/sum(v_init_int_totalR0[which(v_init_int_totalR0$vaccine_intention %in% c(0,1)),]$n)*100

#Vaccine initial intentions per treatment
v_init_int_treatmentR0 <- comboR0 %>%
  group_by(individual_treatment, vaccine_intention) %>%
  summarise(n = n())

vii_pla_1R0 <- v_init_int_treatmentR0$n[v_init_int_treatmentR0$vaccine_intention == 1 & v_init_int_treatmentR0$individual_treatment == "Placebo"][1]/sum(v_init_int_treatmentR0[which(v_init_int_treatmentR0$individual_treatment=="Placebo"&v_init_int_treatmentR0$vaccine_intention %in% c(0,1)),]$n)*100
vii_cdc_1R0 <- v_init_int_treatmentR0$n[v_init_int_treatmentR0$vaccine_intention == 1 & v_init_int_treatmentR0$individual_treatment == "CDC Health"][1]/sum(v_init_int_treatmentR0[which(v_init_int_treatmentR0$individual_treatment=="CDC Health"&v_init_int_treatmentR0$vaccine_intention %in% c(0,1)),]$n)*100
vii_lc_1R0 <- v_init_int_treatmentR0$n[v_init_int_treatmentR0$vaccine_intention == 1 & v_init_int_treatmentR0$individual_treatment == "Low Cash"][1]/sum(v_init_int_treatmentR0[which(v_init_int_treatmentR0$individual_treatment=="Low Cash"&v_init_int_treatmentR0$vaccine_intention %in% c(0,1)),]$n)*100
vii_hc_1R0 <- v_init_int_treatmentR0$n[v_init_int_treatmentR0$vaccine_intention == 1 & v_init_int_treatmentR0$individual_treatment == "High Cash"][1]/sum(v_init_int_treatmentR0[which(v_init_int_treatmentR0$individual_treatment=="High Cash"&v_init_int_treatmentR0$vaccine_intention %in% c(0,1)),]$n)*100

row5 = cbind(vii_pla_1R0, vii_cdc_1R0, vii_lc_1R0, vii_hc_1R0, vii_total_1R0)
row5

sum_tab2 = round(as.data.frame(cbind(t(row3), t(row4), t(row5))),2)
sum_tab2 <- data.frame(lapply(sum_tab2, as.character), stringsAsFactors=FALSE)
colnames(sum_tab2) = c("All Subjects", "Phase 2 & 3 - Compliers", "Phase 2 & 3 - Non-compliers")
rownames(sum_tab2) = c("Placebo", "Health", "Low Cash", "High Cash", "Total Sample")

require(stargazer)
stargazer(sum_tab2,type = "latex", summary = FALSE, out="ExtData_Tables3_PanelB.tex")


# Summary Table 3: 
tab2_working = round(as.data.frame(cbind(t(row3), t(row4), t(row5))),2)

sum_tab3 = round(tab2_working[c(1:4),] - rbind(tab2_working[1,],
                                               tab2_working[1,],
                                               tab2_working[1,], 
                                               tab2_working[1,]),2)
sum_tab3$Difference = round(sum_tab3[,3]-sum_tab3[,2],2)
sum_tab3 <- data.frame(lapply(sum_tab3, as.character), stringsAsFactors=FALSE)
colnames(sum_tab3) = c("All Subjects", "Phase 2 & 3 - Compliers", "Phase 2 & 3 - Non-compliers", "Difference")
rownames(sum_tab3) = c("Placebo", "Health", "Low Cash", "High Cash")

stargazer(sum_tab3,type = "latex", summary = FALSE, out="ExtData_Tables3_PanelC.tex")



